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Abstract 

We present simulation results for a model polymer melt, consisting of 
short, nonentangled chains, in the supercooled state. The analysis focuses on 
the monomer dynamics, which is monitored by the incoherent intermediate 
scattering function. The scattering function is recorded over six decades in 
time and for many different wave-vectors which range from the size of a chain 
to about three times the maximum position of the static structure factor. The 
lowest temperatures studied are slightly above T c , the critical temperature of 
mode-coupling theory (MCT), where T c was determined from a quantitative 
analysis of the (3- and a-relaxations. We find evidence for the space-time fac- 
torization theorem in the /^-relaxation regime, and for the time-temperature 
superposition principle in the a-regime, if the temperature is not too close 
to T c . The wave- vector (q-) dependence of the nonergodicity parameter, of 
the critical amplitude, and the a-relaxation time are in qualitative agreement 
with calculations for hard spheres. For q larger than the maximum of the 
structure factor the a-relaxation time r q already agrees fairly well with the 
asymptotic MCT-prediction r q ~ q~ l l b . The behavior of the relaxation time 
at small q can be rationalized by the validity of the Gaussian approxima- 
tion and the value of the Kohlrausch stretching exponent, as suggested in 
neutron-scattering experiments. 
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I. INTRODUCTION 

During the past decade, numerous experiments and simulations have focused attention 
on the dynamics of supercooled liquids in a temperature region about fifty degrees above the 
calorimetric glass transition |],|2]|||4],|5|] . This interest has been elicited by the development of 
the so-called mode-coupling theory (MCT) ||,§,[7| • Mode-coupling theory predicts that there 
is a critical temperature T c above T g , where the dynamics of the glass former qualitatively 
changes from a liquid-like to a solid-like behavior. The onset of this change is manifested 
by a two-step decay of dynamic correlation functions, which couple to density fluctuations. 
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These correlation functions can be measured in experiments and computer simulations. 
The distinguishing feature of the theory is that it makes universal, system-independent 
predictions about the shape of the time correlation functions, that there are constraining 
relationships between various theoretical quantities, and that key, nonuniversal parameters 
can be expressed in terms of the glass former's static structure. 

Especially the latter point opens the possibility of a direct comparison between theory 
and experiment or simulation. Such comparisons are, however, predicated upon having accu- 
rate data for the temperature dependence of the static structure factor at hand. Therefore, 
they have only been performed for few, simple systems, such as hard-sphere-like colloidal 
particles |J or computer models of soft spheres ||, and binary Lennard- Jones mixtures |jTD| . 

On the other hand, most tests of MCT have concentrated on the universal predictions by 
adjusting the system-specific parameters. This fit procedure must simultaneously optimize 
at least three free parameters. It is rather involved and further complicated by the fact that 
the theoretical predictions are only valid asymptotically in a narrow temperature interval 
around T c , and that the microscopic (vibrational) time scales have to be well separated from 
those of the structural relaxation. These drawbacks have led to a criticism of the significance 
of the experimental evidence for the theory [[□]], but also to extensions of MCT to tackle 
problems, like the interference of vibrational and relaxational time scales [T2|], orientational 
degrees of freedom [13], and corrections to the asymptotic behavior [14.15]. This situation 



suggests that further tests may be beneficial. 

In the present paper we want to discuss the results of a molecular-dynamics (MD) simu- 
lation for a polymer melt and the analysis of the incoherent intermediate scattering function 
by the idealized mode-coupling theory. Comparable applications of the theory to polymer 
melts have already been performed in experiments and simulations. One of the earliest 
applications were neutron scattering experiments of polybutadiene [|16[| . These experiments 
found the predicted square-root behavior for the temperature dependence of the nonergod- 
icity parameter below T c , suggested an in-phase variation of the coherent a-relaxation time 
with the static structure factor, and provided evidence for time-temperature superposition 
above T c . These studies mainly focused on the a-relaxation. A detailed line shape analysis 
in the /^-relaxation regime was not made. Such an analysis was attempted in dielectric re- 
laxation experiments of poly(ethylene terephthalate) fl7j . The experiments are complicated 
by the fact that poly(ethylene terephthalate) has a high tendency to crystallize, and that 
a /3-peak masks the MCT /3-process. Nevertheless, the experiments confirmed some of the 
idealized predictions, but also claimed that there are severe deviations, especially below T c . 
Close to and below T c deviations are expected due to ergodicity restoring processes. Usually, 
the extended MCT [ITS] is then applied to account for these processes. However, a recent 
reanalysis |19[ of the dielectric data indicates that the deviations cannot be explained in 
this way. Instead, a higher order scenario, an ^-singularity, has to be used. In applications 
to other structural glass formers these higher-order singularities are less common, but some 
studies suggest that they could be pertinent to partially crystalline polymers |20| , |2l| . 

On the simulation side, the glassy behavior of polymer melts has been studied by Monte- 
Carlo simulations of the bond-fluctuation lattice model [0,[^,^4[]. These simulations deal 
with short, nonentangled chains, and have been restricted to temperatures above T c . A com- 
parison with MCT indicated that a quantitative description of the /3-relaxation is possible if 
the extended theory is used [BH[. Despite this agreement, there are still some points which 
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could be improved upon. On the one hand, the simulations were performed at constant 
volume, whereas experiments are usually done at constant pressure, and the underlying lat- 
tice structure precludes all phonon contributions to the short-time dynamics. On the other 
hand, the simulation data were not completely equilibrated. They still exhibited very slow 
physical aging processes. These aging processes effectively correspond to additional relax- 
ation channels which are not contained in the idealized MCT. Although it is interesting that 
the extended MCT is able to describe such a situation - especially when taking into account 
the current intensive research on physical aging below T c p6| , |27[| -, the simulations do not 
meet the theoretical premise of thermal equilibrium. 

Both drawbacks are removed by the present MD-simulation. It is done at constant 
pressure and equilibrated on all length scales which we will discuss, i.e., at low temperatures 
the simulation time spent for equilibration exceeded by one order of magnitude the time, 
over which actual dynamical measurements were performed after equilibration. Some details 
of the model, of its static properties, and of the simulation are compiled in the next section. 
The subsequent section is split into two parts. The first describes the analysis of the j3- 
relaxation, and the second that of the a-relaxation. The last section summarizes the main 
results. 

II. MODEL, STATIC RESULTS, AND SIMULATION TECHNIQUE 

This section briefly reviews some characteristics of the model. A detailed description of 
its properties and of the simulation technique can be found in Ref. |28| . 

The simulations were done with linear, nonentangled, monodisperse chains of length 
N = 10. All monomers interact by a truncated and shifted Lennard- Jones (LJ) potential, 
C/ L j(r) = 4e[0/r) 12 - (cr/r) 6 ] + C. The constant C = 0.00775 assures that the potential 
vanishes if r > 2r m i n = 2 x 2 1 / 6 er, where r m j n is the minimum position. Temperature and 
distances are measured in units of e/k-Q and a, respectively, and time is measured in units 
of (mcr 2 /e) 1//2 , where the mass is set to unity. 

In addition, there is a FENE-potential between bonded monomers along the backbone 
of a chain, i.e., E/p(r) = — 15i?Qln[l — (r/R ) 2 } with R = 1.5. With these parameters the 
superposition of the FENE- and LJ-potentials generates a steep effective bond-potential 
with a minimum at about 0.96a 

In our study a melt configuration contained 120 polymers, and ten configurations were 
simulated at each temperature to improve the statistics. The simulation procedure con- 
sisted of two steps. First, the volume was allowed to fluctuate in a simulation at a given 
temperature and pressure (always p = 1 - the influence of pressure is studied in Ref. ) 
to determine the equilibrium density at this thermodynamic state point. Fixing then the 
resulting volume the subsequent runs used the Nose-Hoover thermostat to simulate in the 
canonical ensemble. All dynamic properties calculated are results from these canonical sim- 
ulations. In order to equilibrate the system, each chain was propagated several times over 
the distance of radius of gyration before starting the analysis. This time suffices so that the 
incoherent scattering functions have decayed to zero for all wave-vectors (except for q < 2 
at T = 0.46). 

The particular choice of the bonded and nonbonded potentials has two main consequences 
for the static properties of the model. First, the chains do not become stiffer with decreasing 
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FIG. 1. Temperature dependence of the melt's structure factor. The temperatures span the 
interval from the normal liquid (T = 1) to the supercooled state (T = 0.46) of the melt. The 
lowest temperature is slightly above T c ~ 0.45. In the mode-coupling analysis mostly q = 3, 6.9, 9.5 
are used. The smallest g-value probes the size of a chain (end-to-end distance: R c ~ 12.3, radius 
of gyration: R g ~ 2.09), whereas the larger wave- vectors correspond to intermonomer distances. 

temperature. In the interesting temperature region (T < 0.7) the end-to-end distance, R e , 
and the radius of gyration, R g , are essentially constant: R% = 12.3 ± 0.1, R g = 2.09 ± 0.01. 
Second, the monomer distance 0.96, favored by the bond-potential, is incompatible with 
multiples of r min , and thus with crystalline ordering. This is illustrated by the static structure 
factor of the melt in Fig. [lP. At small wave-vectors, the structure factor is of the order of 
0.01, indicating that the melt has a low compressibility, then it raises and exhibits an 
amorphous halo close to q — 6.9, corresponding to the nearest-neighbor packing of the 
monomers, and finally decreases again and begins oscillating to gradually approach 1 for 
large q. This structure is characteristic of the liquid state and present at all temperatures. 
With decreasing temperature the maxima and minima become sharper, and the position of 
the first maximum slightly shifts to larger g-values because the density of the melt increases. 



1 This figure extends the data of Fig. 6 in Ref. [28|. Note that we chose the conventional nor- 
malization in the present case contrary to Ref. p8[ . In Fig. [I], S(q) tends to the isothermal 
compressibilty divided by that of the ideal gas, i.e., to pTnx, if q — ► 0, whereas it approaches 
(pTnx /total number of monomers) in Ref. [28]. 
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III. IDEALIZED MODE-COUPLING ANALYSIS 



Close to the critical temperature mode- coupling theory predicts a two-step relaxation 
behavior for dynamic correlation functions that couple to density fluctuations, like the in- 
coherent intermediate scattering function, 

W = m £ \ 6XP {[q ' [T " m(t) ~ rm( ° )]) / ' (1) 

m=l ' 

where M stands for the total number of monomers in the melt. The expected two-step 
relaxation gradually develops for T < 0.55 in our model. Therefore the following analysis 
focuses on this temperature region. It is divided into two parts. First, we discuss the 
/3-relaxation regime, and then the final structural decay, the a-relaxation. 



A. /3-relaxation regime 



The /5-relaxation regime is defined as the time window, where the correlator (p s Jt) is close 



to the nonergodicity parameter f q c 
/| c take the following form |14j 



i.e., \<jf (t) - /„ sc | < 1 H|,0,|1|,|T|. The corrections to 



+ h s g G(t) + hf > 



V 2>, 
T J 



(2) 



The first two terms represent the space-time factorization theorem. The spatial dependence 
is contained in f q c and in the critical amplitude h s q , whereas the /3-correlator G(t) carries 
the whole time and temperature dependence. In the idealized MCT it is given by [Bl|l4,15 



G(t) = J\e\g(t/t £ ) with t £ 



t 



o 



\l/2a 



C 



Tr-T 



(3) 



where to is a matching time to the microscopic transient, t e is the /3-time scale, e is the 
separation parameter (C = constant), and g(t) is the temperature independent (3 master 
function. Its shape is determined by the exponent parameter A, which in turn fixes the 
critical exponent a and the von-Schweidler exponent b via A = T(l — a) 2 /T(l — 2a) = 

r(i + 6) 2 /r(i + 26). 

The third term of Eq. (0) is a leading-order correction (of order to h s q G(t) (which is 



of order lei 1 / 2 ) 



It extends the description of the decay from to longer times (i.e. 



for t e <^ t <^ t it is a correction to the von-Schweidler law t ), and depends on temperature 
by the a-time scale r 



to 



e t 



for T>T C 



(4) 



with 7 = l/(2a) + 1/(26). 

When applying these formulas to the simulation data we proceeded in two steps. First, 
we tried to fix the exponent parameter by working with the factorization theorem only. To 
this end, we calculated the /3-scaling function for a specific value of A numerically. The 
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FIG. 2. Incoherent intermediate scattering functions versus time at T = 0.47 (symbols) and at 
T = 0.52 (lines). Three q- values are shown: q = 3 (~ size of the chain), q = 6.9 [~ maximum of 
S(q)], and q = 9.5 first minimum of S(q)]. 

result was inserted in Eq. (Q), and the remaining parameter, f q c , t e , and the total prefactor 
of g(i), h B jT) = /i^l 1 / 2 , were adjusted at a given temperature and g-value. This procedure 
was repeated for different A- and q- values (mainly q = 6.9, 9.5), and at several temperatures 
(mainly T = 0.48, 0.52) to explore which range of A-values yielded fits of comparable quality. 
The best-fit result was: 

A = 0.635 ±0.025 , a = 0.352 ± 0.010 , b = 0.75 ± 0.04 , 7 = 2.09 ± 0.07 . (5) 

In addition, this first step also helps to find the upper bound of the temperature interval, 
where the asymptotic formulas of MCT can be applied. For the model studied this is 
T w 0.52. Therefore, the second step was restricted to T < 0.52. In this step, we fixed 
A = 0.635 and optimized the remaining parameters, /| c , t e , h s q , and h s ^ 2 \T) = h q ( - 2 * ) T~ 2b . 
This procedure was done for T = 0.47, 0.48, 0.49, 0.5, 0.52, and at about 20 different q- values, 
ranging from q = 1 to q — 19, for each temperature. 

The simulation results for T = 0.47 and T = 0.52 are compared in Fig. 0. The figure 
shows the incoherent intermediate scattering function over six decades in time@ for three 
different wave vectors. The smallest g-value (q = 3) probes the length scale of a chain, 
whereas q = 6.9, 9.5 approximately correspond to the maximum and the first minimum of 



. Six decades in time 
therefore correspond to 5 x 10 6 MD steps. 
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Time is measured in units of the MD time step, which is dt = 0.002 [28 




FIG. 3. Incoherent intermediate scattering function versus time at T = 0.47. Three ^-values are 
shown: q = 3 (~ size of the chain), q = 6.9 [~ maximum of S(q)], and q = 9.5 [~ first minimum 
of S(q)]. The dashed lines are the MCT-fit results, using only the factorization theorem, whereas 
the solid lines also include the corrections to the von-Schweidler law. The dashed horizontal lines 
indicate the fit values for f q c at the respective momentum transfers. 

the structure factor S(q) (see Fig. ||). At a given q- value the simulation data almost coincide 
for both temperatures if t < 10 _1 . This short-time regime corresponds to the ballistic motion 
of a monomer, which only shows a weak y/T temperature dependence. Contrary to that, 
the curves strongly separate from one another with increasing time. At T = 0.47, <p s q (t) 
reaches zero about a decade later than at T = 0.52. To achieve a similar growth of the 
structural relaxation time when cooling from higher temperature to T = 0.52, one has to 
take T « 0.7, a temperature that is 35% larger than T = 0.52. Compared to this difference, 
the 10% disparity between T = 0.47 and T = 0.52 is indicative of the approach to T c , where 
the asymptotic formulas, Eqs. (@)-(§|), should hold. 

Figures ^] and |3] show comparisons between the simulation data and Eq. (^) for T = 0.47 
and T = 0.52, respectively. For both temperatures the description of the data by the theory 
starts at about the same time, t ~ 0.6, (except at T = 0.52 and q = 3, where the fits extends 
to t < 10 _1 ). This suggests that the scale to, which MCT introduces as a matching time 
of Eq. (H) to the transient microscopies close to T c (i.e., to ~ to(T c ) = constant), is only 
weakly - if at all - temperature dependent for our polymer model. The same conclusion 
could have also been drawn from Fig. § because the scattering functions, calculated at 
different temperatures, coincide for small times (i.e., for t < 0.3). Whereas the asymptotic 
result, t w constant, appears to hold quite generally in experiments and simulations of 
nonpolymeric systems [Q], a pronounced temperature dependence was found in a Monte 
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FIG. 4. Incoherent intermediate scattering function versus time at T = 0.52. Three ^-values are 
shown: q = 3 (~ size of the chain), q = 6.9 [~ maximum of S(q)], and q = 9.5 [~ first minimum 
of S(q)]. The dashed lines are the MCT-fit results, using only the factorization theorem, whereas 
the solid lines also include the corrections to the von-Schweidler law. The dashed horizontal lines 
indicate the fit values for fi^ at the respective wave-vectors. 

Carlo simulation of the bond- fluctuation model [p5 ]. Compared to the present study, this 



seems to be a mo del- specific rather than a typical feature of glass forming polymers, contrary 
to the conjecture of Ref. p5| . 



For t > 0.6, the idealized theory describes the decay of the correlators over about 1.5 
decades in time at T = 0.52 and over more than 2 decades at T = 0.47. Therefore the 
/^-window expands considerably in this narrow temperature interval. The extension of the 
window can also be inferred from the shift of the time r co , where fq C crosses the simulation 
data. For all g-values it increases from about 1.4 at T = 0.52 to about 10 at T = 0.47. 
The independence of r co on q is an evidence of the factorization theorem which implies 
4> s q ( T co) = fg C for all q, if G(r co ) = 0. As T — > T c + , one expects G(t) = to occur close to the 



/3-relaxation time so that r co oc t e fL5| . For T > T c , t e marks the crossover of the critical to 
the von-Schweidler dynamics, where corrections are not dominant yet. This can be seen in 
Figs. and (|. The idealized fits with and without corrections coincide on the scale r co , but 
deviate at later times. Depending on q, corrections to the von-Schweidler law extend the 
fits by about 0.25 decades at T = 0.52 and by about 0.5 decades at T = 0.47. For q = 3, 
they are negative, but positive for q = 9.5. The change of sign occurs approximately around 
q = 6.9. This behavior qualitatively agrees with theoretical calculations for hard spheres 
14| and recent simulations for linear molecules ||30|| . 

Despite this agreement, we want to mention a problem that we faced when applying 
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FIG. 5. Temperature dependence of the /3-time scale, t £ , and of the total prefactor of the 
/^-scaling function, h s q = /iglej 1 / 2 . Plots of t~ 2a and of (h q ) 2 versus T should extrapolate to zero 
at T c [see Eq. (J|)]. The fit results (dashed line for t £ and solid lines for h q ) can be combined to a 
common estimate: T c = 0.450 ± 0.005. 

the corrections. In the analysis we treated h s q and h q ^ as temperature- and g-dependent 
fit parameters. However, contrary to h s q , the resulting temperature dependence of h q ^ 
turned out to be rather irregular and sometimes even opposite to the theoretical expectation, 
especially at low temperatures. From Eq. (0) one expects the significance of the corrections to 
diminish as T c is approached, whereas Figs. |3] and |] show that they become more important 
with decreasing temperature. This discrepancy could have two reasons. First, h B ( 2 > is about 
2 to 3 orders of magnitude smaller than h s (which itself is of the order 1CT 2 ). Therefore it is 
numerically difficult to extract reliable values from the fits, especially at large q- values, where 
f q c is small. Second, the idealized MCT is only supposed to work in a temperature interval 
which is close, but not too close to T c . In the immediate vicinity of T c , ergodicity restoring 



hopping processes start to compete with the ergodicity breaking cage effect that is 

treated by the idealized MCT. Qualitatively, the main influence of the hopping processes is 
to accelerate the decay of the correlators at late times. Since the uncorrected idealized fits 
for small q lie above the simulation data, the addition of corrections to the von-Schweidler 
law could partly mimic the effect of hopping. This ambiguity can complicate the analysis. 
Presumably, it would have been better to impose the theoretical temperature dependence 
on h s ( 2 \ and to adjust only the variation with q. 

Despite this proviso, the idealized analysis should yield reliable results on the /3-time 
scale around r co , where neither corrections, nor hopping processes are dominant. From the 
fits one can extract the g-dependences of f q c and h s q , and the critical temperature. These 
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FIG. 6. g-dependence of the nonergodicity parameter f q c and of the critical amplitude h q (in 
fact, hgty is shown. The microscopic matching time to cannot be separated in the fit). The solid 
curves are the Gaussian approximations of Eq. (^). The fit interval for f q c was 1 < q < 8. The 
result of this fit was used to estimate /i ms d from the behavior of fit at small q. 



results are shown in Figs. |5] and [(]. Figure |5] plots t~ 2a and {h s q ) 2 versus temperature. The 
theoretically predicted linear behavior [see Eq. @] is nicely confirmed by the data. A linear 
regression yields T c ~ 0.451 for t £ . Similar estimates are obtained for h q : T c ~ 0.451 
(q = 3), T c ~ 0.448 (q = 6.9), and T c ~ 0.448 (q = 9.5). Combining these results we find 
T c = 0.450 ± 0.005, which implies 0.04 < (T — T c )/T c < 0.15 for the reduced distance to the 
critical point in our model. Although (T — T c )/T c = 0.04 corresponds to the upper limit, 
where some of the asymptotic laws begin to be detectable in theoretical calculations for hard 



spheres [|T4],|T5| , similar reduced distances are not unusual in practical applications of MCT 
(see Refs. |U],|32|, for instance). 

Figure |6| shows the g-dependences of f q c and h s q , together with the Gaussian approxima- 
tions, 



J a 



exp[-<? 2 r s 2 a 



and 



h 8 q w h msd q exp[-g r. 



SCJ ' 



(6) 



which are supposed to hold for small q As in other experimental |33| and theoreti- 

cal studies pD| , p4] , the nonergodicity parameter of our model monotonously decreases with 
increasing q, whereas the critical amplitude passes through a maximum between the max- 
imum and the first minimum of S(q). The Gaussian approximations provide reasonable 
descriptions for the initial behavior of both /? c and h q . The description extends farther 
for (q < 10) than for h s (q < 4), which is qualitatively comparable to the theoreti- 
cal results for hard spheres 0,B3]. From the fits one obtains: r sc = 0.095 ± 0.005 and 
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FIG. 7. a-scaling plot of <p q (t) versus t/r q for q = 6.9. Nine different temperatures are shown: 
T = 0.46, 0.47, 0.48, 0.49, 0.5, 0.52, 0.55, 0.6, 0.7 (from left to right). The scaling time T q was defined 
by <pq(T g ) = 0.3 (lower horizontal dashed line). Note that the scaling extends to shorter rescaled 
times with decreasing temperature if T > 0.48. For T = 0.46 (dashed line) this trend is violated 
(see text for further discussion). 



^msd^o = 0.0045 ± 0.0010. Both parameters enter the first two terms of the short-time ex- 
pansion for the a-process of the (monomer) mean-square displacement, which corresponds 
to the von-Schweidler law in reciprocal space [14]. The mean-square displacement will be 



briefly discussed in the next section and more extensively together with related quantities 
in Ref . M . 



B. a-relaxation regime 

The final structural decay from the nonergodicity parameter to zero is called a-relaxation 



regime. For this regime the idealized MCT makes the following predictions |l4] , p!5| , |35[] . First, 
the a-process asymptotically obeys a time-temperature superposition principle, 

f q {t) = 4>l(t/r) , (7) 

i.e., all correlators, measured for one q- value at different temperatures, should collapse onto 
a master curve if time is rescaled by r [see Eq. (Q)]. Second, the short-time expansion of 
Eq. (0) coincides with the long-time behavior of the /5-process, i.e., with the von-Schweidler 
law, (j) s q (t) = fq C — h s q B(t/T) b (B = 0.476 for our model). Third, the leading-order corrections 
to Eq. (0) extend the description beyond f^ c at short rescaled times. They are given by 
64> s q (t/r) = h q (Bi/B)\e\(r/t) b (with Bi = 0.185), and are identical to the corrections of the 
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/^-correlator for t ^> t £ . Forth, in the limit of large q the a master curve is given by a 
Kohlrausch function [p7|1 



with = f* c , Pk = b, and r^(T) = f^r(T), where oc q 1 ' b . On the other hand, it is 
found that Eq. (g) represents a good approximation for the a master function even if q is 
small [|35J. 

Figure shows a test of these predictions for q = 6.9. A scaling time r q was defined by 
< t>q( T g) = 0-3- This is legitimate because any time belonging to the a-regime is expected to ex- 
hibit asymptotically the same temperature dependence, i.e., r q ~ r, due to time-temperature 
superposition. The figure shows that the time-temperature superposition property is borne 
out by the simulation data for T < 0.7. Therefore, it already starts a higher temperatures 
compared to the /3-scaling. Theoretically, such a difference can be rationalized by the fact 
that the corrections to the (3 master function are of order e 1//2 , and thus larger than those 
of Eq. (|7|), which are of order e only. At T = 0.7 the a-scaling is realized for about the 
last 40% of the decay. With decreasing temperature it extends to smaller rescaled times in 
such a way that it progressively adjusts to the von-Schweidler law, and deviations at shorter 
times are quantitatively described by the /3-correlator. 

This qualitative trend persists as long as T > 0.47. However, for T = 0.46 the rescaled 
correlator moves away from the von-Schweidler asymptote instead of approaching it further. 
This is in contrast to the theoretical expectation. Although we could not propagate the chain 
over several times the radius of gyration, we believe that the melt is equilibrated over the 
distances probed by q = 6.9 so that the deviation is presumably not a residual nonequilibrium 
effect. A possible explanation is that the melt at T = 0.46 ((T — T c )/T c = 0.02) is already 
so close to T c that the influence of ergodicity restoring processes starts to be felt. These 
processes change the shape of the correlator in the von-Schweidler regime and prevent the 
strong increase of the a-relaxation time according to Eq. (|J) from continuing |Tj|. It is 
therefore possible that a correlator falls above the idealized a master curve when trying to 
rescale it by a relaxation time which is too small. However, a comparable overshoot was not 
observed in a recent MD-simulation of diatomic molecules ||30|| , although the temperature 
dependence of the relaxation times deviated markedly from the idealized prediction at the 
lowest temperatures studied. Nevertheless, the low temperature curves start crowding in 
the plateau regime and no longer expand along the von-Schweidler asymptote for shorter 
rescaled times. Therefore, it is conceivable that an overshoot could also occur in this model 
at still lower temperatures. 

Figure |] shows that the relaxation time at T = 0.46 is indeed smaller than expected from 
the mode-coupling fit. The figure depicts a plot of t~ 1 ^< and D 1 ^ (D: diffusion coefficient 
of a chain) versus T, using the best estimate for 7 from the /3-analysis (7 = 2.09). In the 
temperature interval, in which the idealized analysis was possible (0.47 < T < 0.52; at 
T = 0.46 the /3-analysis with the parameters quoted above did not yield acceptable results), 
the simulation results lie on straight lines, as predicted by MCT. However, for T = 0.46, 
the relaxation times are too small (no reliable estimate of the diffusion coefficient could be 
obtained after 5 x 10 7 MD-steps). The same observations were also made in Ref. 0] for 
comparable reduced distances to the critical point. Therefore the interval 0.47 < T < 0.52 




(8) 
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FIG. 8. Plot of T q ^ 7 versus temperature for q = 3 (~ size of a chain), q = 6.9 [~ maximum 
of S(q)], and q = 9.5 [~ first minimum of S(q)]. The exponent 7 was taken from the /^-analysis 
(7 = 2.09). The straight lines are fit results for the interval 0.47 <T< 0.52. The intersections of 
the fits with the zero-line are estimates for T c : T c = 0.446 (q = 3),T c = 0.449 (q = 6.9), T c = 0.452 
(q = 3). The inset shows the same plot for the diffusion constant. Here the fit result for the critical 
temperature is significantly smaller, i.e., T c = 0.438. 




was used to obtain further estimates for T c . For r q this yields T c m 0.446 for q = 3, T c m 0.449 
for q = 6.9, T c w 0.452 for q = 9.5, and for D, T c w 0.438. Whereas the variation of T c for 
r q , albeit decreasing systematically with decreasing q, is compatible, within the error bars, 
with T c = 0.450 from the /3-anaiysis, the result from D is too small. Barring the problem 
that the diffusion coefficient is hard to determine accurately at low temperatures (i.e., for 
T < 0.47), the high temperature behavior (i.e., 0.48 <T< 0.52) suggests that the disparity 
in the T c -estimates is significant. Such a finding is not uncommon. Similar observations 
were made in MD-simulations of a water model [38| and of a binary Lennard- Jones mixture 
Whereas the difference in T r in the water simulation was small and could be attributed 



to numerical uncertainties, a much larger reduction was found for the binary mixture. On 
the other hand, an unconstrained fit yielded a critical temperature which was compatible 
with the estimates from coherent and incoherent scattering, but 7 was in turn significantly 
smaller than that from the /^-analysis. The same result was also obtained in MD-simulations 
of diatomic molecules p0|j39[1 and for our model [f29| . 

In experiments and simulations it is generally found that the Kohlrausch function rep- 
resents a very good fit formula for the major part of the a-relaxation. Only at early times 
deviations (can) become visible. Mode-coupling theory interpretes these deviations as a 
signature of the von-Schweidler law, since it concatenates the (3- with the a-process, and 
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FIG. 9. Kohlrausch stretching exponent /3k [see Eq. as a function of g for three different 
choices of the fit interval 0.02 < (jf q {t)/f s q c < X for X = 0.4, 0.7 and 0.85. The thick horizontal line 
indicates the von-Schweidler exponent b = 0.75. In addition, three different q- values are highlighted 
by vertical dashed lines: q = 2i:/R g 4.35 (R g = radius of gyration), q = 6.9 [« maximum of 
S(q)}, and q = 9.5 first minimum of S(q)]. 



P q 7^ b in general 

arises that the late /3-process interferes with the fit 



Therefore, when applying the Kohlrausch function, the problem 

Mj, and additionally, that the fit pa- 



rameters sensitively depend on the size of the fit interval (for a discussion of this problem 
see Refs. ||4~I] , |4"2fl , for instance). To lessen this problem two procedures are often used. On 
the one hand, one can fit the late-time decay of the a master function by taking into ac- 
count that the short-time bound of the fit interval should not overlap too much with the 
von-Schweidler regime. This procedure has been used in Refs. p0| , p2| , for instance. On the 
other hand, one can work with so-called a-/3-fits, in which Eq. (0) (without corrections) and 
Eq. (0) are superimposed. Except for hard-sphere colloidal systems, where it is possible 
to calculate the master function numerically Eq. (^) is usually replaced by Eq. 

The superposition requires = fq C , and the von-Schweidler law to be subtracted from the 
/3-correlator because it is already approximated by the Kohlrausch function (for an 

alternative procedure see plf). 

A variant of the latter approach was applied in this simulation. We only worked with 
the Kohlrausch function, but fixed the amplitude by = f^ c . Then we tried to optimize 
the stretching exponent /3k (the time scale is already given by <j> s q {r^) = /| c e _1 ). The results 
are still very sensitive to the choice of the fit interval. An example is given in Fig. [5| It 
shows the g-dependence of the stretching exponent Pk- The percentages specify the portion 
of the decay from /| c to (about) 0, which is included in the fit (i.e., "40%" means the range 
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FIG. 10. 4>g(t) and Kohlrausch functions at T = 0.48. Solid lines are the simulation data for 
nine different g-values: q = 1,2, 3, 4, 6, 8, 10, 12, 14 (from right to left). For q > 2, the Kohlrausch 
parameters are given by: /? = f q c , /3k = b = 0.75, and from <fr s q (T^) = f q c e~ l . These choices 
provide a good description of the ce-decay, as long as q > 2. For q < 2, deviations gradually 
develop. The stretching becomes more pronounced, leading to the following parameters at q = 1: 
= f q c , 0k = 0.656, and T q = 24595. In addition, two vertical lines are shown. The thin dashed 
line is t e (T = 0.48) = 4.933, whereas the thick dashed line indicates the lower time value (t = 500), 
above which the monomer mean-square displacement exhibits a i°' 65 -behavior (see Fig. |IT| ). 



0.02 < cj) s (t)/ f q c < 0.4). At large q the nonergodicity parameter is so small that statistical 
noise of the data considerably influences the results. The /^-values for the various fit 
intervals then strongly splay out, leading to a difference of about 15% at q — 19. With 
decreasing q the curves approach one another and fluctuate around (3k ~ b = 0.75, except 
for the smallest g-values, q = 1,2, where /3k lies between 0.65 and 0.7. 

This behavior suggests that it should be possible to find a fit interval for q > 2, where 
the Kohlrausch function with = and (3k = b provides a good description of the a- 
process. Figure [1(] shows a comparison between 4> s q (t) and such Kohlrausch functions for a 
selection of g-values, ranging from q = 1 to q — 14, at T = 0.48. One can clearly see that 
the expectation /3k = b is well borne out if q > 3, but that smaller stretching exponents are 
required for q < 3. 

In addition, the data exhibit a small bump at about t = 3 for all wave vectors. This bump 
gradually develops at low temperatures (T < 0.52), but is hard to see in the presentation of 
Figs, fj or [| A similar feature was also observed in other simulations of fragile glass formers, 
for instance, for orthoterphenyl or for a binary Lennard- Jones mixture J34], and is much 
more pronounced for strong glasses, like silica |4"6|| . Usually, the bump is attributed to a 
finite size effect. Due to nonlinear coupling with other wave vectors a sound wave, which 
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propagates through the simulation box, leaves it on one side, and immediately reenters the 
box on the other side because of periodic boundary conditions, could generate a disturbance 
at a time t »s L/c, where L is the linear dimension of the box and c is the sound velocity. 
Using L = 10.5, and estimating the sound velocity as c ~ 7,1 we find t ~ 1.5, which is close 
to the position of the bump. 

Two comments with respect to the Kohlrausch fits have to be made. First, the equality 
/3k = b starts to work at g-values which are about a factor of 2 smaller than the maximum 
position of S(q). Such a behavior has neither been observed in experiments [f44| , f45"f nor in 
simulations |32j of nonpolymeric glass formers. Theoretically, one expects (3k = b to become 
valid if q — > oo, and numerical calculations for hard spheres indicate that this limit is ap- 
proximately realized only for q- values which are 5-6 times larger than the maximum position 
of S(q) fl3~5| , |3"7| . Therefore the present finding seems to be a special property of the model 
under consideration, although neutron-scattering experiments for polybutadiene [|l6l , [47| and 
other polyolefmes f|S],[|U| also suggest /3k = const (and perhaps ~ b for polybutadiene |T6] ) 
for q- values that lie around the maximum of S(q) (typically, 27r/i? g -C0.2<g<5A ). 

The second comment concerns the drop of (3k for small q. Usually, one expects (3k to 
increase monotonously to 1 as q — > because the system becomes freely diffusive on this 



length scale. This behavior is well borne out in theoretical calculations for hard spheres 
and simulations of simple liquids ||34|| . We believe that the difference observed here is related 
to the polymeric character of our model because q < 2 probes the length scale of a chain. 
For these q- values the Gaussian approximation, 



(j> s Jt) = exp 



-~? 2 0o(*) 
6 



(9) 



represents a (very) good description of the data (see Fig. [11]) so that one can interprete the 
decay of <p s q (t) in terms of the monomer mean-square displacement^ go(t). At early times 
(t < 0.03) the mean-square displacement is purely ballistic, i.e., go(t) = 3Tt 2 = 1.44t 2 
(T = 0.48), before it crosses over to a plateau regime at about t ~ 0.3. At this time, Eq. @ 
also ceases to be a good approximation for (fi s q (t), if q is too large (i.e., for q > 5). The 
leveling-off of go reflects the temporary localization of the monomers in their cages. It is the 
counterpart of the /3-relaxation in real space, and can be described by |14| 



9o(t) ~ 6r s 2 c - 6h 



msd 



g(t/t £ 
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(10) 



where C a and Cb are constants, which result from the corrections to the critical decay and to 
the von-Schweidler law in the limit q — > 0, respectively. For hard spheres, one expects C a < 
and C b < |14l,[15|. Figure p4]b includes a fit of Eq. (|I~0"D to go(t). Please note that rf c = 0.009 



3 The sound velocity is estimated from c = {c v jcy pnx) 1 ^ 2 by taking into account S(q — > 0) = TpKT 
and the thermodynamic relation between the specific heats, c p and cy, and the thermal expansion 
coefficient. 

4 All monomers are included in go{t), as it was the case in the calculation of 4> q (t). For a further 



discussion of this and other mean-square displacements see Ref. [36| 
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FIG. 11. Panel (a) shows a comparison between <pg(t) (symbols) and the Gaussian approximation 
Eq. (||) (solid lines) for five different (/-values at T = 0.48. The (/-values smaller than 3 probe the 
size of a chain (see Fig. |l]). For q = 1 (< 2ir/R e ) the agreement is quantitative. Panel (b) shows 
the mean-square displacement, go(t), of all monomers (solid line) at T = 0.48. The initial slope 
is g (t) = 3Tt 2 = 1.44t 2 (dashed line). The dotted line is a fit by Eq. (0) with 6r s 2 c = 0.054 
(dashed horizontal line). In addition, a power law, go(t) ~ t x ° with xq = 0.65, is indicated by 
another dashed line. To determine xq, only the simulation data larger than R 2 = 2.09 were used. 
Nevertheless, the power law extends by about one decade to smaller times. The dash-dotted line 
shows the expected long time behavior go(t) = 6Dt. The diffusion coefficient was determined from 
the mean-square displacement of a chain, which reaches the diffusive limit earlier than go(t) |36]| . 
In both panels three vertical dashed lines are shown. From left to right, the first (t = 0.3) indicates 
the approximate time, when the Gaussian approximation stops working, if q > 5, the second is the 
/3-time scale t e (T = 0.48) = 4.933, and the last marks the onset of the t ' 65 power law. 
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(r sc ~ 0.095 corresponds to the Lindemann criterion of melting), h ms atQ = 0.0045, t e = 4.933, 
B = 0.476, and g(t) are taken from the /3-analysis. The only adjustable parameters are C a tg 
and CbtQ. The fit yields the reasonable values C a tg ~ —0.3 and C^q ~ —0.25. On the other 
hand, it is also possible to fit g (t) by the leading-order ansatz g (t) — 6r 2 c + At b , where r 2 c 
and A have to be adjusted |50| (see also Ref. |32| for another application to simple liquids). 
Although such a fit extends the description of Eq. (|ToD at late times by approximately one 
decade, it requires r sc « 0.087 and A « 3.68 x 10~ 3 . These values are considerably different 
from those of the /9-analysis. Therefore, the first approach should be preferred. 

Equation ( [TP] ) describes the simulation data from t ~ 0.5 to about t ~ 100, i.e., up to the 
initial relaxation of a monomer out of its cage. During the j3- and the early a-relaxation a 
monomer hardly feels the bonds to its neighbors along the chain. The mean-square displace- 
ment is of the order 10~ 2 — 10 _1 of the monomer diameter. Therefore the monomer behaves 
like a particle in a simple liquid. For larger displacements chain connectivity starts to influ- 
ence the monomer dynamics, and finally becomes dominant. Whereas the von-Schweidler 
behavior crosses over to free diffusion in simple liquids [j32]| , a further sub diffusive regime, 
9o(t) ~ t x ° (xq ~ 0.65), intervenes for our model if t > 500. Such subdiffusive displace- 



ments are well known in polymer simulations of nonentangled chains [SI], and usually in- 
terpreted as a crossover between a t 1 / 2 -behavior - the prediction of the Rouse model |52 
for 1 <C g (t) <C -R 2 - and free diffusion. A comparison of Figs. and [TT] shows that the 
subdiffusive displacement occurs in the same time regime, where the stretched exponential 
with /3k = 0.656 fits <jf q {t) for q = 1. Therefore, we suggest that the drop of (3k is a polymer 
specific effect, related to chain connectivity, which becomes visible if q probes the length 
scales between the monomer diameter (= 1 in our units) and the end-to-end distance R e . 
On the other hand, if q <C 2%/R e , the monomer displacement has to become diffusive, and 
so /3k should approach 1, as in simple liquids. Therefore we expect /3k to exhibit a minimum 
for our model at small q. 

Finally, Fig. |I2| shows the g-dependence of the a-relaxation time r q . Here, r q was defined 
by 4> q {r q ) = f q c /2 instead of by because rf is larger than the simulated times if q < 2. 
However, for g-values, where both times are available, we checked that they qualitatively 
yield the same result. In addition, r q was divided by t\ +a l b to eliminate the critical depen- 
dence on T c [see Eq. (U)]. Asymptotically, the resulting quantity, f q = T q t [see Eq. (§)], 
should be temperature independent. The figure shows that the division does not completely 
remove the temperature dependence. For q = 1 the relaxation times at T = 0.47 and 
T = 0.52 are about a factor of 2 different, with f q (T = 0.47) < f q (T = 0.52), whereas both 
f q s first approach, and then cross one another around the first peak of the static structure 
factor (q ~ 6.9), before splaying out again. At q — 19 there is again about a factor of 2 
between T = 0.47 and T = 0.52, but now f q (T = 0.47) > f q (T = 0.52). However, the 
absolute values of f q are of the order lO^ 1 for q = 19, and thus 1-2 orders of magnitude 
smaller than tl +a ^ b , while they are comparable to tl +a ^ b at q = 1. Therefore the difference 
at small q seems to be more important. It could again indicate that the relaxation times 
on the largest length scales do not increase as rapidly as expected from idealized MCT, an 



observation, which we have already made for the diffusion coefficient (see Fig. || and Ref. [29 
for further discussion of this point). On the other hand, it is not clear how significant the 
found deviations are at all, since the dominant temperature dependence has been taken into 
account by tl +a ^ b , and an additional smooth temperature dependence is expected, if T is 
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FIG. 12. q-dependence of the a-relaxation time T q , defined by 4> q (r q ) = f q c /2. The figure shows 

T q t aj/b = T q /t\ +a ^ b [see Eq. (|j)] for three different temperatures from the interval, where the 
/3-analysis could be done. Thick solid lines indicate different asymptotic behavior. The initial 
decrease of f q t a ^ b is steeper than q~ 2 , the expected dependence for free diffusion, but compat- 
ible with q~ 2 IPK = ^-3.05 ^ K _ 9555) j see ^|iT|)] . For larger q, f q t a ^ b crosses over to the 
mode-coupling result, Eq. ([l2]). The filled points show (f^ c /Bhg) 1 ^ , whereas the solid line is the 
asymptotic limit, ~ q~ 1 ^ b , for q — » oo. 



not very close to T c (see Fig. 7 of Ref. 15 , for instance) 



The g-dependence of the a-relaxation time can be roughly divided into two regions. For 
q < 2 Figs. |TD| and [IT] showed that the Gaussian approximation, Eq. (|9[), and the Kohlrausch 
function with /3k = 0.656 describe the a-relaxation of 4> s q (t) well. In the a-regime one can 
therefore approximately equate Eqs. and ([]). Furthermore, if one assumes a power law 
for the g-dependence of r q , one obtains 



q 



-2/Pk 



q 



-3.05 



for q < 2 



(11] 



Figure |T2| shows that this estimate reasonably agrees with the initial behavior of r q . Equa- 
tion ([11]) has been suggested already some time ago in neutron-scattering experiments of 

,2, 



glass- forming polymers JIE|,18|], if q go(t) is sufficiently small to warrant the Gaussian ap- 
proximation (see also Ref. [47| for a critical discussion of this issue). If q > 5, f q is indicative 



of another power-law behavior that is compatible with the mode-coupling prediction |35| , |37 

fsc i 1 / b 

J, ?Z^ g -l/6 ) (12) 



Bh\ 



which is expected to hold, if /3k = b. 
experiments of orthoterphenyl |33|| . 



A similar observation was made in neutron-scattering 
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IV. CONCLUSIONS 



This paper presents results of a molecular-dynamics simulations for a simple model of a 
glassy, nonentangled polymer melt. The discussion focuses on the monomer dynamics above 
T c , as monitored by the incoherent intermediate scattering function (jf q (t) and the mean- 
square displacement. These functions do not distinguish between bonded and nonbonded 
monomers so that it is a priori not clear how the polymeric character influences 4> s q (t) in the 
supercooled state. From the analysis the following picture emerges: 

1. At very early times the motion of a monomer is purely ballistic. This corresponds to 
the phonon contribution in laboratory glass formers. 

2. At later times, the ballistic motion slows down because the monomer feels the confine- 
ment imposed by its nearest neighbors. The monomer becomes almost localized, and 
the period of this localization extends over about two decades in time in the narrow 
temperature interval 0.47 < T < 0.52. "Localization" in this context means that the 
mean-square displacement remains close to r sc « 0.095, i.e., to a displacement that is 
about 10% of the monomer diameter. This is the regime of the /^-relaxation, which 
deals with the approach towards and leaving of the plateau value r sc (see Fig. |TT[b). 
In this sense, r sc can be interpreted as a dynamic measure of the size of the cage. Mo- 
tions on this length scale are so small that they are hardly affected by the polymeric 
character, i.e., by the nature of the glass former. This is presumably the reason, why 
mode-coupling theory, a theory developed for simple liquids, can give a reasonable 
explanation for the dynamics of structurally much more complicated systems in the 
/^-regime. 

3. Concerning the comparison with the idealized MCT in the /3-regime, the simulation 
provides evidence for the space-time factorization theorem with temperature indepen- 
dent f* c , h q , and A, as long as 0.47 < T < 0.52. In this temperature region, the 
/3-relaxation time exhibits the predicted power law, yielding a critical temperature of 
T c ~ 0.45, and the g-dependences of f q c and h q are in qualitative agreement with 
calculations for hard spheres. For T < 0.47, deviations from the idealized behavior are 
observed, which can be interpreted in terms of ergodicity restoring hopping processes. 
We have not attempted an extended MCT-analysis because data for only one tem- 
perature below T = 0.47 had been simulated. However, we can roughly estimate an 
upper bound of the hopping parameter 5 for our model. Using 5t = (t /t £ ) 1+2a |L 



t e (T = 0.47) = 8.75, and assuming t « 0.3 (due to Fig. [IT]), we obtain St « 3 x 10" 3 , 
which would be larger than experimental values at similar distances from the critical 
point fnj||2 . 



4. When the mean-square displacement becomes larger than r sc , the monomer begins to 
leave its cage, and the temporarily frozen structure "melts" (Lindemann criterion). 
The initial stage of this "melting", i.e., of the a-process, can be described by the 
von-Schweidler law. But for larger times, polymer specific properties start dominat- 
ing the dynamics. For instance, the mean-square displacement exhibits a subdiffusive 
behavior between the von-Schweidler law and free diffusion due to chain connectivity. 
The exponent of the corresponding power law is typical of the short chains studied. 
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It also determines the Kohlrausch stretching exponent for small wave-vectors, where 
the Gaussian approximation holds. Qualitatively, this difference between a- and (3- 
processes is expected by mode-coupling theory because the space-time factorization of 
the /3-regime, which implies the same dynamics on all length scale, is no longer valid 
in the late a-regime. However, the theory still predicts that the diffusion coefficient 
and the a-relaxation time r q should exhibit the same temperature dependence. In 
this respect, we find deviations between theory and simulation. The diffusion coeffi- 
cient does not decrease as quickly as r q at the maximum of the structure factor with 
decreasing temperature. Somehow the melt stays more mobile on large length scales 
than on local ones. On the other hand, we find evidence for the time-temperature 
superposition principle for the a-process, if T > 0.47. 

5. The wave- vector dependence of the a- relaxation time consists of two regimes. For 
small q, where the Gaussian approximation holds, we find a g~ 2 / /3K -behavior, as in 
some neutron-scattering experiments ||I6|J48|| , which crosses over to a g -1 / fe -behavior 



for q > 5. The latter power law is a MCT-prediction for g-values much larger than the 
maximum of S(q). Why this asymptotic behavior can already be observed for rather 
small q in the present model, is not clear. 

In summary, the dynamics of our model in the supercooled state can be understood as an 
interplay of mode-coupling and polymer specific effects. The "caging" of a monomer by its 
neighbors leads to a temporary trapping of the monomers and to a concomittant slowing 
down of the structural relaxation, as in simple liquids. This cage effect is dominant, as long 
as the monomer displacement is small (i.e., much smaller than the diameter of a monomer). 
The polymeric character of the model only determines the nonuniversal parameters of MCT. 
However, if a monomer gradually leaves its cage, chain connectivity becomes more and more 
influential. For the present model of short nonentangled chains, it leads to a subdiffusive 
Rouse-like displacement. The subdiffusive behavior interferes with the late-/3/early-a dy- 
namics and limits the von-Schweidler regime to a range, which is much smaller than in 
simple Lennard- Jones liquids, before free diffusion sets in. 
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